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ABSTRACT 

In this paper we design and develop several filtering strategies for the analysis of data 
generated by a resonant bar Gravitational Wave (GW) antenna, with the goal to assess 
the presence (or absence) in them of long duration monochromatic GW signals, as well 
as their eventual amplitude and frequency, within the sensitivity band of the detector. 
Such signals are most likely generated in the fast rotation of slightly asymmetric 
spinning stars. We shall develop the practical procedures, together with the study 
of their statistical properties, which will provide us with useful information on each 
technique's performance. The selection of candidate events will then be established 
according to threshold-crossing probabilities, based on the Neyman-Pearson criterion. 
In particular, it will be shown that our approach, based on phase estimation, presents 
better signal-to-noise ratio than the most common one of pure spectral analysis. 

Key words: Gravitational Waves - resonant detectors - monochromatic signals - 
data filtering 



1 INTRODUCTION 

It is generally believed that the most intense Gravitational Waves (GW) arriving in the Earth from remote sources in the 
Universe correspond to very short duration (~1 millisecond) bursts, generated in the explosion of a supernova (Thorne 1987), 



or in gamma-ray bursters (Roland et al. 1994). Since their very first origins, cylindrical bar GW antennae have been applied 
to the detection of this sort of events ( Weber 1969|), and the mo re modern cryogenic bars have been used for this purpose, 



too, with considerably enhanced sensitivities (Astone ot al. 1991; Astone et al. 1993; Hamilton et al. 1994): the long decay 



times of the oscillations of the bar make it well suited for the measurement of impulsive, short duration signals (Gibbons & 



Hawking 1971 Astone, Bonifazi fc Pallottino 1990) 



It so happens however that some cylindrical GW antennae have been in continuous operation regime for many consecutive 
months, even years. This is the case for example with the Explorer detector, owned and operated by the ROG group at Roma 
(Italy) and installed within CERN premises in Geneva (Switzerland) (Astone et al. 1993). Long term operation naturally 



provides the appropriate background for a search of monochromatic signals in the detector data, as requisite long integration 
times become available. 

Monochromatic signals are most probably generated by the rotation of asymmetric stars, such as a pulsar or a neutron 
star. The intensity of the GWs strongly depends on the amount of asymmetry of the source, and this is in turn dependent on 



its equation of state (Bonazzola & Gourgoulhon 1996). Reasonably optimistic upper bounds on typical star parameters give 
an extremely weak signal estimation of ft~10~ (Thorne 1987), which must be seen against a noisy background. Clearly, 



long integration times are required to reveal this kind of signal. 

A systematic search for it must face a practical difficulty which derives from the fact that the signal is received in the 
antenna Doppler-shifted due to the daily and yearly motions of the Earth — in addition to possible internal motions within th e 
source if it is e.g. in a binary system. Fourier analysis of long stretches of data results in high frequency resolution ( Kay 1990 ), 
thence in signal spread across several spectrum bins if it is Doppler shifted. This can naturally cause significant reduction 
in post-filter signal-to-noise ratio. The problem is easily overcome if the source position in the sky is known (or assumed) 
ahead of time by means of suitable corrections based on ephemeris calculations. Analyses of this type exist in the literature: 
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traces of a pulsar in the centre of the supernova SN1987A were sought in four days of data generated by the 30 metre 
Garching interferometer in March 1989 (Niebauer et al. 1993), and Frasca and La Posta studied almost four years of data 
generated by the room temperature bar detector GEOGRAV in search for periodic signals from the Galactic Centre and the 
Large Magellanic Cloud (Frasca & La Posta 1991). More recently, Mauceli (Mauceli 1997) has looked for monochromatic GW 
signals coming from the region of Tuc J^l and from the Galactic Centre in three months of data generated by the cryogenic 
detector ALLEGRO at Luoisiana State University. 

A different strategy must of course be used for an all-sky search. The philosophy of the procedure put forward by Frasca 
and La Posta (Frasca & La Posta 1991) consists in the construction of a large bank of spectra, taken over shorter stretches of 
data such that the frequency resolution in each individual spectrum be sufficiently low that daily Doppler shifted signals fit 
in a single spectral bin. Suitable comparison and averaging are thereafter applied to the spectra in order to draw conclusions 
about the intensity and/or bounds on signals. Astone et al. (Astone et al. 1997; Astone 1998), have looked at one year (1991) 
of data taken by the above referenced Explorer detector to also perform an all-sky search for monochromatic sources of GWs. 
Their method is based upon local maxima identification in a bank of spectra, followed by close up analyses of frequency peaks 
looking for evidence of Doppler shift patterns across the duration of the entire data stretch. 

In this paper we design and develop algorithms for the analysis of data generated by a resonant bar detector of GWs, in 
search for monochromatic signals within the system's sensitive frequency band. We are also interested in an all-sky search, 
but adopt a different point of view. Rather than scanning a bank of spectra, we propose to use a matched filter technique to 
estimate both frequency and phase of candidate signals, then set a threshold, using the Neyman-Pearson criterion, to select 
those events which have a given probability of crossing it as a consequence of pure random noise fluctuations. We have tested 
our methods in simulations with real Explorer detector data from 1991, and seen that they perform very satisfactorily. We 
plan to apply them to the massive processing of long stretches of data from the same antenna in a future paper, in order to 
provide complementary analyses to the procedures and methods already reported in (Astone et al. 1997) and (Astone 1998). 

The article is structured as follows: in section 2 we present a few technical generalities and set the basic notation 
conventions. Section 3 is devoted to a detailed study of a situation in which the signal has a frequency exactly equal to one of 
those in the discrete Fourier spectrum of the data (Lobo & Montero 1997); this corresponds to an idealised situation whose 
consideration is methodologically useful, as it allows us to determine the signal's phase, and to investigate the statistical 
properties of the filter output; it also characterizes the main guidelines for the more realistic study in subsequent sections. In 
section 4 the method is illustrated with an artificially added signal to real detector data, which includes the estimation of the 
noise spectral density in the presence of such signal. In section 5 we address the real case, in which the signal frequency no 
longer exactly matches any of the discrete samples, so that it leaks across neighbouring spectrum bins (Lobo & Montero 1998), 
and also assess the statistical properties of the filter output (Montero 1998). Finally, in section 6 we apply the method again 
to real data with an external control signal added, and show that it works satisfactorily. The paper closes with a summary of 
conclusions and future prospects. 



2 LINEAR DATA FILTERING 

We begin with a review of some fundamental concepts of linear data processing, fixing also the basic notation which we will 
be using throughout this article. 

In the general case, let u(n) (n = 0, . . . , N — 1) be the discrete set of samples which constitute our experimental data. A 
linear filter consists in a discrete set of numbers g(n; fii) depending on several parameters, fa, which acts on the experimental 
data as follows: 

N-l 

producing what we shall call the filter output. It is usually assumed that u(n) is the sum of two different contributions: on 
the one hand the signal, x(n), whose presence we want to assess, represented by a deterministic function, and on the other 
hand the noise r(n), a stochastic process: 

u(n) = x(n) + r(n). (2) 

For any choice of parameters, it is appropiate to ask for the filter response both to the signal, yx, and to the noise, yr, the 
latter being a stochastic process, too. The ratio of the mean square values of these quantities is termed in the literature the 
output signal-to-noise ratio (SNR), 



and it is a measurement of the performance of the filter g(n; fa). The theory of the matched filter (Helstrom 1968; Papoulis 
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1984) precisely determines, up to a global constant, the functional form which this must have, for given signal and noise, in 
order to maximize p. 

In our case we shall assume that the noise r(n) can be adequately modelled by a zero-mean Gaussian and station- 
ary stochastic process, whereas the signal x(n) will be the response of the cryogenic resonant detector Explorer to a pure 



monochromatic GW (Pallottino fc Pizzella 1984; Montero 1997) 



x(n) = Ao cos(2nfonT + ipo). 



(4) 



Here, Ao is the product of the amplitude by a conversion factor which defines the detector sensitivity at the frequency of the 
gravitational radiation, f g . This differs from fo by a constant shift ( Frasca et al. 1992] ), 



fim =fg-fo = 900.0267781 Hz, 



(5) 



introduced by the data-acquisition system of the antenna, with the purpose to sample the antenna's full bandwidth 27.5088 
Hz. The matched filter for such a signal is then functionally equal to the latter (Helstrom 1968), 



g(n) = Bcos(2-k f nT + ip ), 
where B is an arbitrary constant. 



(6) 



3 NON-LEAKING SIGNALS EMBEDDED IN KNOWN SPECTRUM NOISE 
3.1 The signal 

As pointed out in the introduction, we want to develop in this article a general method whose operativity does not depend on 
the existence of prior information about the source. So, in principle, the value of the frequency fo of the signal we can detect 
must be within the interval 

0</o<^, (7) 

where 1/T is our Nyquist frequency. Obviously no search strategy can afford the endless analysis of all the frequencies in 
that window, so we shall be forced to select a finite set of frequencies to scan. Nevertheless, the very functional form of the 
filter shows us that we shall perform discrete Fourier transforms (DFTs) in its implementation, and this defines the set of 
frequencies which will be searched in actual practice: 

fT = jf ke{0,...,N-l}. (8) 

Moreover, for practical reasons, all the DFTs will be numerically computed using the fast Fourier transform (FFT) algorithm, 
a very optimized procedure which naturally computes at once all spectral components, with the only restriction that the 
number of samples be an exact power of 2. 

In this section, we shall assume that the signal is well matched by the spectral template. By this we mean that the 
frequency of the signal is in fact one of those in equation (^), so that all the signal is in one single bin of the FFT, with no 
leakage to the neighbouring ones. More precisely, we shall be assuming that 

foT=^, (9) 

where ko is one of 1, . . . , f - 1, though we do not know which. We shall disregard the study of any ko bigger than N/2 
because it would be redundant since they represent nothing but negative frequencies. The value fco = is also disregarded 
since, among other considerations, it represents no wave at all, but a constant signal. 

Summing up, the target of the present analysis will be to assess the presence of a signal 

x(n) = A cos(2nk n/N + tp ) (10) 

in the experimental data series u(n), using a matched filter 

g(n; k, ip) = B cos(2n kn/N + ip). (11) 

depending on the two unknown paramenters, k and ip, which we shall eventually estimate. Besides the advantageous property 
of the absence of frequency leakage in the filter output of such signals, equation (j^|) shows that x(n) is a periodic function 
over the entire processed period, because x(iV) is equal to x(0). In fact, this relationship holds for any sample, 

x(n + iV) =x(n), (12) 

and it will be a crucial aspect for the developments which we shall introduce below. 
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3.2 The filter performace and the role of B 

Let us compute the two quantities y x and < y£ >, in order to evaluate the actual goodness of the filter. y x is different from 
zero only if a signal is really present and the value of the parameter k matches ko, i.e. the filtered signal does not leak across 
different frequency bins, 



y-x(k,<p) = 



cos (ip — ipo) Skk - 



For < yr(k, <p) > we have the more complex formula, 

2 JV_1 / l\ 

<yl{k,p) >=^- R(n) (l - jj-) cos(2nkn/N) - B 

n=-(JV-l) ^ ' 

where we have introduced the autocorrelation function of the noise, 
R(n) = (r(n + n)r(ra')}- 



2 cos(27rfc/7V - 2ip) 



sm(2nk/N) 



^ R(n) sm(2nkn/N), 



(13) 



(14) 



(15) 



If we assumeu that R(n) ~ for n ~ iV, it is clear that the second term shall become negligible in front of the first, and 
hence, 



<y r (fc,<^) >«F 



2 jV S(fc; iV) 
2T 



with 



JV-l 



5(fc; N) =T ^ J?(n) f 1 - U. j cos(27rfcn/7V). 

n=-(JV-l) ^ ' 

It can be easily shown that the quantity we have just defined is the mean value of a periodogramme, 

JV-l 



-i2irkn/N 



(16) 



(17) 



(18) 



a well-known way for estimating the power spectral density of the noise at that particular frequency, based on the Wiener- 
Khinchine theorem (Kay 1990). 

Putting together expressions ( |l3| ) and (|l^) we finally get for the SNR, 



A NT 2/ \ r _ 2/ N r 



with 
Po = 



AqNT 
2S{k;N)' 



(19) 



(20) 



the maximum value for p we may achieve with the present filter. 

The SNR is obviously independent of the constant 8, so we can freely set it at our will in order to provide y with some 
advantageous property. Our particular choice is, 



B(k) 



2T 



NS(k;N)' 



(21) 



the factor that makes < yr(k, <p) > equal to one. The statistical properties of the noise and the linearity of the filter guarantee 
that y is still Gaussian. Then its probability distribution will be completely settled once we know its mean < y >, and its 



variance, cry which in our case coincides with < y r >, 

a y =< y > - < y > =< yr >= i- 



(22) 



So, on the one hand, we have forced ay to take the same value regardless of the particular scanned frequency, and on 
the other hand, the mean of y(k, ip), 



< y(k, ip) >= yfp^cos(ip - ipo) S k k , 



(23) 



* This, in fact, is a constraint on the number of samples that we want to filter. We must set N large enough to sufficiently exceed the 
correlation time. 
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shall be zero when either no signal is present in the data or, if there is a signal!!, for any value of k other than fco. In this way, 
we have designed a bank of filters whose outputs corresponding to pure noise are statistically equivalent, and consequently 
can be directly compared. 



3.3 Data splitting and averages 

So far, we have implicitly assumed that N represents the total amount of stored information we have access to or, in other 
words, that we can analyze the whole data series in a single filter pass. This is, in many senses, a rather optimistic assumption. 
First of all, since we want to process several months of experimental data, it should not be surprising that the available (or 
even existing) computing facilities could not afford such a calculation. Moreover, the output of any experimental device, like 
Explorer, will not be uniform in quality along all the data acquisition time, and the stationarity of the noise is not preserved 
over too long periods of time. Then it could be worse to mix bad data (those, for instance, with a high level of noise) with 
good data in a single analysis, than simply to veto the stretch that we find unacceptable. But the gaps that we may introduce 
in rejecting samples of the experimental set are not the unique discontinuities that we shall find in the time series, because in 
such a long term operation of a detector it is not unlikely that the system suffers sporadic stops. Also the properties of the 
physical signal could be not so stable to be satisfactorily fit by our models along extended periods of time. 

We shall thus consider that each series of length N is just one among a set of, say, M consecutiveH blocks (Astone et al 



1997; Astone et al. 1997). The reasons for the choice of the particular values of N and M must not necessarily coincide in 
general. In particular, it is possible that there exist several of those sequences of iV x M data, eventually disconnected, which 
must be then processed separately. 

So, we shall attach a new label a to each quantity in order to be able to specify which of the M blocks of N data we refer 

to: 



y a (fc, tp) = J fc- — J2 u " ( n ) cos(27rfc7i/iV + tp) (o = 0,...,Af-l), (24) 

n — 

a is actually a shorthand which simplifies the notation, 



u a (n) = u(n + aN). (25) 

It is obvious, however, that computing y a (k, tp) for each a will not change the individual values of p, as shown by equation 
(jil^). Our final goal should be then to combine them in a suitable way which allows us to make the final SNR as high as 
possible. The definition ( |25| ) is, in this sense, very revealing because, when combined with ([lo]) and ([l2|), shows the most 
important feature of a non-leaking monochromatic wave: the signal x a in fact does not depend on a, 

x a (n) = Ao cos(2nkon/N + tpo) = x(n). (26) 

We consequently see that, if the signal is present, each of these blocks contains an identical replica of the same stretch of 
sinusoid in them. This motivates us to define a new random variable z(k,tp), 

M—l 

< k >^ = ji X) y °(*'*0' ( 27 ) 

whose mean value does not differ from that in equation (^3|) , but whose variance is reduced as a consequence of this averaging 
operation. Before we calculate explicitly this quantity and the value of the new associated SNR, we are going to focus on the 
problem of choosing the right value of the phase parameter of the filter, tp. 

The conceptually simplest method is to compute z(k, tp) for a lot of different values of that parameter, then select the 
best, tp, i.e., that which gives a larger output after the filtering procedure. 

Nevertheless, we do not need to go into such computationally long process, for the optimum value tp can be analytically 
determined as follows. According to its definition, we may write down z(k, tp) as 



{k >^ = M\ NKk-N) £ l«{Mfc)} cob p + ^Mk)} amid, (28) 



1 / 2T 

where 



t In fact, wo will put together both cases, since we shall take the criterion of typifying through po = any frequency that contains no 
signal, whether or not there is a signal at some other value of fc. 

t The relaxation of this condition, allowing for the existence of missing whole blocks, introduces minor changes in all the following 
discussion. So the derived expressions can be easily adapted to this case. 
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u a (k) = j2 u ^ e ~ 27V ' kn/N ( 29 ) 

71 = 

is the DFT of u a (n). So, we define the (p imposing a local-maximum condition on z(k, ip): 
dz(k, ip) 



d(fi 



= (30) 



We thus finde 

M-l 



^3{u«(fc)} 

tan(£) = j£ . (31) 

^5R{u a (fc)} 

Here it is useful to introduce the two random variables, IR(fe) and I(fc), 



, M-l 

1 / 2T \ - 



M V ^^(fc; JV) ^ 

a— 



K{u Q (fc)}, (32) 



K*) - ^V^i&v)^ 9{G4fc)} ' (33) 

CK — 

because they can easily be combined to yield z(k, <p) 7 



z(k) = z(k, (p) = vfiHlf, (34) 

where we have defined z(fc) as the best output of the filter at a given frequency, extending the notation used with tp. 

The statistical properties of the actual filter output z(k) will strongly differ from those of z(k, <p), since the new random 
variable is the fruit of a non-linear filtering process. For instance, its mean is no longer equal to zero, even if there is no signal 
in the experimental data, as we shall see. This is the reason why we did not undertake a very detailed study of z(k, tp) in the 
first place. 

3.4 Probability distribution of z(k) 

We shall build the probability density p(z) starting from p(K) and Those two auxiliary random variables are Gaussian 
by construction, and are statistically independent, as they are the real and imaginary parts of the Fourier transform of 
a stationary stochastic process. Then we only have to know the respective means and variances in order to complete the 
information that will fully settle their probability densities. The mean values of K and I can be readily found from their 
definitions and equation (p3|): 

<R(fco)> = Vpo"cos^ , (35) 
<I(A;q)> = ^/po"sin<^o (36) 
while for the variances it can be found that both quantities are approximately equal, and 

o-r w o-, 2 » — , (37) 

because the correlation time is much shorter than the duration of the individual series, which is in essence equivalent to state 
that, in spite of their consecutiveness, they are mutually uncorrelated. 

The hypotheses we have made lead us to the following expressions for the probability density of R and I: 



p(R) = \ /ilL e --(«-VPoco S¥ > , and ^ 
p(I) = y^ e -¥C-VPo-«.) a . (39) 



We can ensure that ip defined in this way leads to a maximum and not a minimum of function z(k, tp) because its second derivative, 



d%,z(k,tp = tp) = — z(fe), is never positive. 
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Hence p(z) is given, after an integral is evaluated, by 
p(z) = A/ze-* (z2+Po) /q(MzVpo), (40) 



Jo(l/)=i/ e-" -"^ (41) 



o 



where we have used one of the integral representations of the modified Bessel function of order zero ( Olver 197C ), Io(y) 

In the absence of signal, i.e., when po = 0, equation (|40|) reduces to 
p(z) = Mze~^ z2 , (42) 

an expression which explicitly displays the statistical equivalence of all the frequencies which contain no signal, the property 
we want to achieve when we set the value of the constant B. Moreover, in this case, z is Rayleigh distributed or, in other 



words, z follows a \ distribution with two degrees of freedom (Papoulis 1990) 



3.5 Mean and variance of p(z). A new SNR. 

We are now interested in the mean and variance of z. These correspond to the first two moments of the probability distribution 
p(z). It appears that a closed analytic expression can be found for the moments of any order, so we consider it here for 
completeness. 

The m-th moment is defined by 



<z m >= / Mz m+1 e-- (z +Po) / (Mz^), (43) 
Jo 

a calculation that becomes straightforward if one uses the relationship ( fcradshteyn fc Ryzhik 198C| ), 



Ly(-z) = jr^yy J x*e-^I {2^Tz)dx, (44) 
where L y (z) is the Laguerre function of order y, assuming the normalization condition, 

L y (0) = 1, (45) 
and r(y) is Euler's gamma function. The expectation value of the m-th power of z is then, 

<*™>=GI)M? +i M-f4 <«> 

The most relevant moments for our purpose are, as has been said, the first and the second. The mean, 
as we announced, will be different from zero even when po vanishes, due to the property ( ^j|) of Laguerre functions, 
Nevertheless, the asymptotic behaviour of the Laguerre function 

U-'^W+T)*> (49) 

also shows that < z >, when a signal is present, approaches the maximum mean value which according to ( jj^ ) the random 
variables y a {k,tp) can possibly reach, 

< z >^ ^/po~. (50) 



It is worth noting that equation ( |49[ ) ensures that the last expression holds not only when po 3> 1 but when p> 1, where 
we have introduced 

a quantity that plays the role of the new SNR. The same conclusion can be obtained after the study of the explicit expressions 
for z| : 

Zx = po (52) 
and < z r >, 
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< 4 >= A. (53) 

It is relevant to point that, in this case, SNR linearly increases with the total number of filtered samples, N x M: 

_ AqNMT .... 

P =TS(^N)' (54) 

whereas the more classical procedure of averaging the square of the moduli of the DFTs leads to 

_ a 2 nVmt 

Po = — po = AS{k . N) , (55) 



what means that with our strategy for signals whose frequency is one of the FFT samples, we enhance by a factor vM the 
value of p. 



4 NON-LEAKING SIGNALS EMBEDDED IN UNKNOWN SPECTRUM NOISE 
4.1 Replacing S(k; N) 

It is almost redundant to say that the operative method we have just developed requires knowledge of the power spectral 
density of the noise. The aim of this section is the effective substitution of S(k; N) in the definition of the constant B by a 
suitable estimate of this quantity obtained from the same data series. 
Let us begin with a rearrangement of expression ([l8|), 

S(k;N) = ^ <[(K{u tt (fc)}-K{x(fc)}) 2 + (9{u Q (fc)}-9{x(fe)}) 2 ]). (56) 

There is one procedure in this formula that is certainly beyond our control: we cannot perform the statistical average. 
Our particular choice will be the substitution of that operation by a sum over the the entire rank of values of a, because 
S(k;N), in spite of the formal aspect of (|56|), is independent of the block label. The same applies upon replacement of x(fe) 
(obviously also an unknown quantity) since x(fc) =< u(k) a >, 

M-l 

*(*) — » jjj? £ Q «w- ( 57 ) 

a=0 

Summing up, the random variable we shall use in order to estimate S(k; N) is S(fc; N), 

\ 2 ' 

(58) 



M-l 



M—1 \ 2 / M-l 



K{Q a (fc)} - i K{u a ,(k)} + Q{u a (k)} ^ u -'( fc )} 

V a'=0 / \ a'=0 J 



where we have divided by M—1 and not by M because this way we get an unbiased estimator, i.e., 

< S(k-N) >= S(k;N). (59) 

Now we can replace the unknown spectral density by E>(k; N)in any preceding expression, thus obtaining a new filter 
output Z(fc) which, unlike z(fc), we are able to compute directly from the raw experimental data u a (n). With an analogous 
procedure to the one already explained we obtain p(Z) and all its related quantities, including the corresponding SNR. Instead 
of starting from scratch, we shall calculate the probability density of Z in two steps, using previous results. 

Let us introduce the auxiliary random variable W(fc) 

W,*) = ||||, (60, 
which allows us to define Z(fc) in a simple way, 

Z=-^=, (61) 



thanks to the fact that all the terms containing S(k;N) in both random variables mutually cancel out. Since W and z are 
statistically independent, and p(z) was given in the last section, we have thus reduced the problem to obtaining p(W) and 
performing a final integration. 



The probability density of W can be found in most reference books on Probability (Papoulis 1990), because it is the 
arithmetic mean of the squares of 2M — 2 zero-mean independent Gaussian variables with unit variances. So (2M — 2) W 
follows a x 2 distribution with precisely 2M — 2 degrees of freedom, 

P(W) = ^-/^p W^-V^- 1 ^. (62) 
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4.2 The actual filter output Z and its distribution 



The ratio in ( |6l[ ) which we have used for defining Z is a familiar one in elementary statistics, and is known to follow a Student's 
t-distribution if z is a zero-mean Gaussian variable. Nevertheless, the expression for p(z) is far from a normal density function, 
as we have shown in (KtOj), what compels us to perform an explicit calculation, leading to the result 



p(Z) = 



MZ 



[1 + 



,z 2 



■ exp 



Mp 



M y2 



La 



M 2 Z 2 p 



(63) 



4M - 4 + 2A/Z 2 

where once again we have made use of the formula (^). When no signal is present in the data at one particular frequency, 
expression ( ps] ) reduces to 

MZ 



p(Z) = 



[1 



T z 2 



(64) 



2M-2 

As a matter of fact, Z 2 in this case follows a Fisher's ^-distribution with 2 and 2M — 2 degrees of freedom, because it is 
the ratio of two \ 2 random variables with those degrees of freedom, respectively. 

In order to compute the moments of the density function of Z, the simplest approach is not to use the final expression of 
p(Z) but an intermediate formula 

_ \Af-l 

p(Z) 



MZ 



(M-iy 



f M - 1 e- wlM - 1+ ^ z] I (MZy^oW)dW, 



r(M- 1) 

that will avoid the problem of the integration of Laguerre functions with negative arguments. We thus find 



<z m >= 



(M-iy 



-+i 



r (m 



r(M) 



+ 1 la 



A I 



-po 



(65) 



(66) 



where we have chosen a layout that emphasizes the resemblance with the result corresponding to < z m >. The term inside 
the square brackets approaches unity when M 3> m, and then we recover the formula (ff6|). It is especially interesting to note 
that, in particular, the newly defined SNR remains unchanged. Let us split the second order moment of p(Z) 

1 2 /. M 



< z > 



M - 2 M 



(l + T Po 



in two terms, namely < Z r > and < Zx >, 
< Z 2 >=< Zr > + < Zx > . 



(67) 



(68) 



When no signal is present, po = 0, the value of < Z > is merely due to the response of the noise to the filtering procedure, 



< Z r >= 



M - 1 2 



M-2M' 

so we will accordingly assign to the signal the rest of the outcome, 
1 



< /x >" M^2 
and therefore, 



Po, 



<Z X > 
< Z r > 



M 

T Po. 



(69) 



(70) 



(71) 



4.3 A practical example 

This section is devoted to show the result of such procedure when applied on a small stretch of data taken by the Explorer 
in August 3rd of 1991, and successive days. The starting date was randomly selected since the final purpose of the present 
analysis is not extracting conclusions on the presence of GWs but on the practical performance of the method itself. We have 
thus externally introduced a sinusoidal signal with the required absence of frequency leakage in order to check the ability of 
the method for revealing it. The signal, corresponding to a GW with an amplitude of ho = 10~ 23 , was placed at about 921.4 



Hz, near the detector's plus resonance (Astone et al. 1993). For this particular date the level of noise in the detector was such 
that the SNR for this signal was po ~ 1/3 for a number of filtered samples of iV = 131072, a little less than forty minutes. The 
signal was therefore completely buried in the noise. The specific value iV = 131072 may seem arbitrary in this context, but it 



has a physical reason: it ensures that no Doppler shift can be observed in the individual blocks of N samples ( Astone et al 



1997 ). Once we set iV we can pin down the precise frequency bin which contains the control signal: ko — 50918. 

Let us see what happens when we process six hours of data (M = 9). Looking at Figure [l].a it is by all means impossible to 
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(a) 



a. 1.5 



(c) 



0.5 



M=9 




910 915 920 925 
Frequency (Hz) 



900 905 910 915 920 925 
Frequency (Hz) 



(b) 



2.5 



B 1-5 



(d) 



0.5 



900 



900 905 910 915 920 925 
Frequency (Hz) 



905 



910 915 920 
Frequency (Hz) 



925 



Figure 1. Output of the analysis procedure when the filtered data stretches extend over: (a) six hours, (b) twelve hours, (c) a day and 
(d) two days. Every time series begins at the same instant of August 3rd of 1991. The arrows point to the signal in each of the last three 
cases. Even though they are progressively more prominent, only in (d) is the signal the highest peak. In (c), for example, it is hidden 
between two taller noise fluctuations. 




decide whether the signal is really present or not: SNR has only risen to a value near unity from the original 1/3 with such few 
blocks processed. By increasing M, i.e., processing longer stretches of data, SNR grows and the signal becomes progressively 
more distinct, as we see in Figures [jjb-d, corresponding to half a day, one day and two days of filtered data, respectively. We 
must stress at this point that the theoretical prediction of equation (^T|) that energy SNR grows linearly with the number M 
of processed blocks is very accurately observed in real practice, as we have numerically verified with the plotted data. The 
improvement by a factor of y/~M relative to more standard procedures — see equation — is thus firmly established not 
only in theory but also in actual practice. —. 

With the output of the filtering process for the values of k other than ko we can computed the distribution of Z, because 
when no signal is present it does not depend on k. The case M = 9 is again very interesting because it offers us the possibility 
of comparing the experimental distribution with p(z) and p(Z), thus checking that Z really follows the second and not the 

^ In fact, we do not need to remove ko (or any other presumed signal) necessarily before we compute the output distribution, because 
at most it is only one point in 2 16 . 
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Table 1. Value of the threshold Ao corresponding to different choices of Qo for each of the considered block number M. 



So 


M = 9 


M = 18 


M = 36 


M=72 


icr 3 


1.56 


0.97 


0.65 


0.45 


icr 5 


2.39 


1.35 


0.87 


0.59 



first — see Figure ^.a. For higher values of M the two probability densities converge and become almost indistinguishable 
from one another (Figure |^.b) . In both instances, however, the coincidence between theoretical and experimental distribution 
is remarkable. 

Moreover, the explicit form of p(Z) is very useful, not only in order to compare it with the experimental one, but to fix 
a threshold Ao, which shall help us in the task of deciding whether a given crossing has statistical significance not. We will 
calculate the error of the first kind, or false-alarm probability, Qo, as a function of Ao, 



Qo 



P (Z)dZ : 



; A 2 



(72) 



and then we shall set the upper bound depending on the number of false alarms (i.e. mistakes) we can afford, using the well 



known Neyman-Pearson criterion (Helstrom 
Ao = 



\ 



2M 



1 



M 



- 1 



• So 



(73) 



Table [l] shows how Ao varies both with respect to the value of Qo and the number of processed data blocks, M. The 
signal's height, in units of the plots in Figure is 0.473, 0.684 and 0.595 for M — 18, 36 and 72, respectively. It is therefore 
above threshold if M = 72 with false alarm probability Qo = 10~ 5 , and both if M — 36 and M — 72 with false alarm probability 
Qo = 10~ 3 . Even so it cannot be clearly told from other random fluctuations, as we see in the Figure. This is because it is 
very weak, of course. We shall come back to the discussion of the significance of these thresholds below. 



5 LEAKING SIGNALS EMBEDDED IN KNOWN SPECTRUM NOISE 
5.1 A leaking signal 

We are going to start this section considering the effects that the presence of a general frequency signal in the data may 
produce in the results we have exposed in the preceding sections. In particular, we shall study the new statistical properties 
of the random variables y a (k, tp), that determine the characteristics of R and I, and consequently of z. 
So, in the following we shall relax the condition (^) 

k + £o 



foT = 



N 



(74) 



by introducing the real parameter eo 
1 1 

-2^ 0< 2' 

whose effective consequence is that x a (n) is no longer independent of a, which appears in the form of an accumulative phase 
shift whenever eo is different from zero: 



(75) 



x a (n) = Ao cos(27r[fc + e ]n/N + 2nae + <p ). 



(76) 



While it is true that the frequency remainder also lowers the maximum filter output, due to the spectral leakage of the 
signal, 

sin(7re ) 
' N sm(ne Q /N) 

it is nevertheless the block dependence that damages the filtering procedure, because it is responsible of the sinc-like behaviour 
of the mean values of R(feo) an d I(fco), when considered as functions of M: 



■ COS I 7T£o 



+ ifo ~ <fi + 2naeo 



(77) 



> 



I(fco) > 



sin(neoM) 
NM sin(7reo/A0 

sin(neoM) 



^7re 



NM sin(7re /iV) 
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M ~N. 



M 



1 

N 



+ <Po 



+ Vo 



(78) 
(79) 
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This fact turns the search strategy not so robust as desired, in the sense that for a given eo other than zero, there always 
exists a value Mo (~ 1/eo) for M, above which the SNR decreases noticeably. The frequency band where the analysis method 
works efficiently thus decreases with the number of averaged blocks, a very undesirable feature. We now investigate how this 
problem can be addressed. 



5.2 Phase- varying filter 



As already stated, since the origin of the problem is a carried over phase, we shall solve it introducing a new block-dependent 
filter with one more parameter, e: 



2T 



NS(k;N) 



cos(2irkn/N + ip + 2irae). 



(80) 



with the purpose to compensate for the phase shift. Starting from this new g a {n; k, ip, e), we can calculate the value of each 
y a (k,ip,e), 



y a (k,ip,e) = 



2T 



-|u a (fc)| cos(4' cl (k) — if — 2nae), 



NS(k; N) 

where the following notation has been used: 
Ua(fc) = \u a (k)\e^ k '> . 



(81) 



(82) 



Through a definition formally identical to (|27|), we shall establish z(k, tp, e). Once again it is possible to obtain <p using a 
local-maximum condition, like that in vRu), 



dz(k, if, 



dtp 



u Q (fc)| sin($ Q (fc) — 2nae) 



= 0, 



tan(y;) = 



(83) 



^ |u a (Jfe)| cos(<f> a (fc) — 2-Kae) 

Q=0 

The value of (p leads now to the following expression for z(fe, <p~, e), 
1 / 2T ~ 

a relationship that involves a new quantity, u (k; e), which formally is also a DFT, 

M-l 

u(fc;e) = >J u a (k)e~ 27rla \ 

a=0 



(84) 



(85) 



The template for e, just like in the case of the frequency grid, will be dictated by the convenience of the use ofthe FFT 
algorithm in the computation of u (fc; e). We shall therefore estimate eo within the following discrete rank of values!] of e: 



= Jf ge{0,...,M'-l}, 



(86) 



where, in principle, M' , by the way an exact power of 2, must not necessarily coincide with M. M' has to be greater than M 
if we do not want to waste available information, but on the other hand, it seems that larger and larger values of M' should 
produce an endless increment of precision in the estimation of eo- Nevertheless, as we will show in the next section, M' should 
be kept as low as possible. 

Equati on (|35|), t ogether with the span condition bears a considerable formal resemblance with the so called zoom 

transform (Yip 1976; Hung 1981; de Wild et al. 1987), in which the twiddle factor is now identically equal to one. It must 
however be stressed that, in the present context, ours is an interpolation formula rather than a frequency resolution algorithm. 
Our formalism can thus be naturally extended with the purpose to refine the spectral resolution, simply taking M' =M, i.e., 
data sets of N x M' points, where all possible values of q must be considered, for a selected choice of frequency bins, k. 



II As usual, any value for q equal or larger than M' /2 is related with a negative e. 
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5.3 Statistical properties of z 

The mechanism for finding q is then very simple. We must compute all the z(k, <p, q/M') and select that q which gives the 
largest output, just defining z(k) = z(k, ip, q/M'), i.e., 

z(k) — max {z(k,<p,q/M')}. (87) 

<je{o m'-i} 

When no signal is present at some frequency, it can be proved that again, 

p(z) = Mze"* z2 , (88) 

no matter the choices for k a nd q. So, in th at case, the probability density of z(k), since it is the maximum of M' equally 
distributed random variables (Papoulis 1990), is given by 



t M 7^ 

p(z) = MM ze~~ 



A/'-l 

(89) 



We have not been able to compute < z > other than in the form of an alternating finite series, which is almost useless for 
obtaining generic conclusions about it. Instead of the mean value of z we shall compute its most probable value. The function 
p(z) reaches its maximum (when M' 3> \[e) for 



-|lnM', (90) 

a quantity that decreases as M increases, and when M' decreases. This shows the convenience of setting M' as the first exact 
power of 2 greater than M. With the second-order moment of the distribution we have in principle a similar problem, although 
in this case the alternating series can be transformed into a non alternating one, 

2 2v^( M'\ (-1) + 2 1 

i=i v 7 i=i 
For large values of M' we can approximate the result using the definition of Euler's 7 constant, 

<z 2 >^ A (i n A /' + 7). (92) 

Since the output of the procedure in absence of noise is 
sin(Treo) 

ZX = VPO Ar . , 7at\ : (93) 

N sm(7reo /N) 
the new SNR shall be (if M' > 1) 

/ sin(Treo) V M . . 

P ~ P ° \Nsm(ne /N) J 2 (In Af + 7) ' 1 ' 

The standard procedure of averaging the square of the moduli of the DFTs gives for a general leaking signal the following 
SNR: 

i(ne ) \ 2 \/M 



^ = "°Usin(^o/JV) 1 ~' (95) 



what represents an improvement with the present method of 

p_ _ y/M 

po _ In M' + 7 ' 

a ratio larger than one when (^) holds, which increases with the number of blocks. 



(96) 



6 LEAKING SIGNALS EMBEDDED IN UNKNOWN SPECTRUM NOISE 
6.1 A spectral estimator 

The task of replacing S(k; M) in the filter definition is much more complex than in the non-leaking case. The natural starting 
point is the value of u (k; q/M') for q other than q, but this presents two main problems. First of all, the fact that we are only 
able to choose q out of a discrete set leaves open the unpleasant possibility that the value of the signal frequency lies just in 
the middle of a bin. This means that we shall ignore the precise way in which the signal energy will be distributed among 
the different q's, and thus the relative magnitude of u (k;q/M') when compared with the rest of outputs. The scenario can 
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be even worse, because we have no guarantees that q/M' corresponds to the best approach to the actual eo: the noise can 
mislead us into an inaccurate value of the signal frequency, leaving thus the true u (fe; eo) among the discarded ones. 

Instead of constructing an estimator for S(k; M) in the hope that none of the previously stated possibilities really takes 
place, what could lead us again to a filtering procedure too sensitive to the signal's peculiar properties, we are going to choose 
a democratic estimate S(k;N): perhaps it will not be so accurate as it could, but it will not show appreciable differences in 
its performance depending on the actual frequency of the GW. 

We define S(k; N) through an expression that closely resembles that in equation (681), 



S{k;N) 



T 1 

NM M> - 2 



E 

95^9 



5R{u (k; q/M')} - -jj—- J2 ^{u (*; q/M')} 



9'#9 



+ U{u (k; q/M')} - j^-j ]T ^{u (k; q'/M')} 



9V9 



(97) 



where the u a {k) have been substituted by u (k; q/M'), and the sums do not contain the term where the signal is supposed to 
be located, u(k;q/M'). 

Once more we replace S(k; N) by S(fc; N) in the definition of z, in order to get a new random variable Z which we can 
compute using only experimental data. This quantity inherits two characteristic traits from the way we estimate the spectral 
density of the noise. 

On the one hand, if the signal is large enough the filter output may have a saturation limit, which will depend upon the 
particular values of some parameters, such as M or M' . This means that Z will not surpass a certain threshold, even if the 
amplitude of the signal increases indefinitely. The reason for that behaviour can be found in the fact that when the signal is 
much more intense than the noise, what becomes really difficult to assess is not the presence of the former but the properties 
of the latter. So in these cases S(k; N) will overestimate S(k;N). It must however be stressed that the whole effect results 
in a change in the value of po, and thus does not actually set an upper bound in the SNR p. On the other hand, when no 
signal is present at a particular frequency, S(fc; iV) will underestimate S(k; N), since we do not use q when computing it, and 
u (k; q/M') has the biggest modulus among all the u (k; q/M'). 

To compute the overall probability density of Z(fc) is very difficult, both in a single step (even if no signal is present) 
or in two steps, like in Section 4.1, since the auxiliary random variable W(fc) defined like in (^) is no longer independent of 
z(k). We shall thus obtain it by a different method. We introduce a tentative p(Z) inspired by the limiting probability density 
p(z), because when M is large enough both must coincide: 



p(Z) = MM' wie- 



ld' - 



(98) 



The new parameter w, which condenses all the differences between p(z) and p(Z), measures in some suitable sense the 
bias of the spectral estimator S(fc; iV), i.e., 

This point of view is somewhat vain since we are not in a position to compute < S(k; N) > theoretically, as just stated. 
The value of w, however, can be estimated from the filter output itself, using for instance the relationship (|9l|), replacing z 
by y/wZ, and the statistical mean by an average over the filter output: 



1=1 



N/2-1 

m ^2 z2 ( fc ) 



(100) 



The procedure we shall use to compute (100) also reminds us that Z(fc) can be obtained from the experimental time 
series, irrespective of its probability density. We set the value of the w quantity a posteriori, once the filtering process has 
finished. In fact, the functional form of p(Z) has its very origin in the comparative study of the actual distribution Z(fc) 
and the probability density of reference, p(z), for different values for the free parameters. We pursue these ideas in the next 
section. 



6.2 A single filtering process 
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(a) (b) 




900 905 910 915 920 925 0.3 0.4 0.5 0.6 0.7 0.8 0.9 

Frequency (Hz) Filter output 

Figure 3. Here we find (a) the output of the filtering procedure, Z(fc), for N = f31072, M = 36 and M' = 64, what represents about 
one day of data, and (b) its experimental distribution compared with the theoretically expected, p(z), and corrected one, p(Z). 



(a) (b) 




900 905 910 915 920 925 0.25 0.3 0.35 0.4 0.45 0.5 0.55 0.6 

Frequency (Hz) Filter output 

Figure 4. Here we find (a) the output of the filtering procedure, Z(fc), for TV = 131072, M = 72 and M' = 128, what represents about 
two days of data, and (b) its experimental distribution compared with p(z) and p(Z): in this case it is hard to distinguish from the one 
another. 



Let us begin analysing a stream of about one day of the Explorer data with the layout used in Section 4.3, i.e. N= 131072 



and M—36 M' =64. In fact, we will choose exactly the same time series, starting on August 3rd of 1991, in order to be able 
to compare the non-leaking and the leaking methods. So we get the results shown in Figure ^.a, where we have also externally 
introduced a signal, with similar characteristics with respect to the previous example. Its amplitude is now slightly bigger, 
/io = 3xlO~ 23 , i.e., SNR is about 1, because the new method is slightly less sensitive than the non-leaking one, as we have 
shown. And, of course, the signal spreads across different frequency bins: fen = 50918 and eo =0.1. 

A plot of the distribution of Z(fc) is displayed in Figure ^.b, where it is contrasted with p(z ) in the theoretical 

probability density, and with p(Z) once w was computed following the prescription shown in (10C), the corrected one. The 
agreement of the latter with the experimental probability density is again remarkable. 

It is convenient to point out that when the number of processed blocks increases, the value of w rapidly approaches 
one, thus becoming an irrelevant parameter. As a matter of fact, Figure ^, where we present the same procedure with the 
same data, but extending the processed time to two days, shows that p(z) is then a sufficiently accurate expression for the 
probability density of the actual filter output. 

In the present case we can also compute the error of the first kind in terms of the threshold An: 



Qo 



M' 



and to invert this relationship, 



Ao = 



-^(l-tl-fio]*). 



(101) 



(102) 



By way of example, the threshold levels in units of the graphs in Figures ^| and ^ are 0.96 and 0.68, respectively, for a 
false alarm probability Qo = 10~ . In either case the signal is clearly above these thresholds, as it has respective heights of 
1.12 and 0.83. From here we can rather accurately determine eo, too: we find q = 6 when M'=64 (i.e., eo.estimated = 0.094), 
and q= 13 when M' = 128 (i.e., eo.estimated = 0.101) for a real value of eo = 0.1. 
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7 OUTLOOK 

The long term operation of cryogenic GW detectors opens the possibility of looking for long duration GW signals in the data 
generated by them, since signal-to-noise ratios are enhanced by the availability of long integration times. Monochromatic, as 
well as stochastic signals belong in this category, though the latter require data from two or more independent antennae. In 
this paper we have addressed the problem of the design of suitable algorithms to single out possible monochromatic signals, 
coming from any direction in the sky, in the background of the noisy data produced by a cylindrical bar. 

This is not such a simple problem, due to a variety of reasons: computers are not arbitrarily powerful, detector duty 
cycles are not 100%, data quality is not homogeneous, Doppler shifts distort monochromaticity, etc. Although some of the 
characteristics of the data are peculiar to the detector system, there are a number of procedures which should be quite generally 
usable. This paper is concerned with the problem to set up filtering algorithms which enable the selection of candidate signals 
on the basis of threshold crossing criteria of the filtered data. To this end we have considered banks of filters which have a 
sinusoidal form but which enable signal phase estimation. We have determined the pdf 's at filter output and thence consistently 
established the probability of crossing a given level. We have also designed suitable estimators of the noise spectral density 
which take into consideration the possibility that the signal be in one or more of the FFT frequency bins for the complete 
duration of the analysed data. 

These procedures have been checked by means of simulations on top of real data generated by the Explorer detector of 
the Roma group in 1991, and they work quite well: the theoretical predictions on the improvement in SNR as longer stretches 
of data are processed remarkably coincide with those observed in the real data analysis. They thus appear promising, and 
we would expect them to be useful to analyse other resonant detector data, whether cylinders or the future projected large 
spherical antennae, and also the large scale interferometers currently under construction, with suitable adaptive modifications. 

As we have seen with the simulations in this paper, the Explorer detector comfortably sees amplitudes of 10 -23 , and 



therefore one can expect it to be sensitive to signals a few times 10 -24 — see also (Astone et al. 1997). These are still rather 
high values for the expected amplitudes of signals of astrophysical origin, so threshold crossing criteria must be established 
having this in mind: if we set a very low false alarm tolerance then chances of missing the real event will grow, whereas if we 
set it rather low, the real signal will be treated on the same footing as random noise. 

Real signals however have characteristic Doppler patterns, whose correct assessment should provide a more sound selection 
of candidate signals. In a forthcoming article we plan to apply the methodology here developped to a systematic analysis of 
the above mentioned Explorer detector data over a long period of time, and to extend it to also include Doppler shift effects in 
the data. The very algorithms of section |^ will be a powerful tool in this respect, if implemented under the perspective of the 



zoom transform (Yip 1976), to enhance the frequency resolution within the spectral neighbourhood of pre-selected candidate 
lines. 



ACKNOWLEDGMENTS 

The authors are grateful to the members of ROG group at Rome for giving them access to the data of the antenna Explorer, 
as well as the detector's transfer function details and procedures related to specific issues of the data analysis. We also thank 
Guido Pizzella for his critical reading of the manuscript and valuable suggestions, and for his continued hospitality. Bernard 
Schutz's comments, pointing us to the zoom transform algorithm, have been very helpful and are greatly appreciated, too. 
MM wishes to thank the DGR of the Generalitat de Catalunya for financial support, and JAL the Spanish Ministry Education 
for a grant, PB96-0384. We also acknowledge support from the Institut d'Estudis Catalans. 



REFERENCES 

Astone P., Bonifazi P., Pallottino G.V., 1990, Rev. Sci. lustrum., 61, 3899 
Astone P. et al., 1991, EPL, 16, 231 
Astone P. et al., 1993, PRD, 47, 2 

Astone P., Frasca S., Pallottino G.V., Pizzella G., 1997, in Bassan M., Ferrari V., Francaviglia M., Fucito F., Modena I., eds, Proceedings 
of the XII Italian Conference on General Relativity and Gravitation Physics, World Scientific, Singapore 

Astone P., Frasca S., Pallottino G.V., Pizzella G., 1997, in Ciufolini I., Fidecaro F., eds, Proceedings of the VIRGO International 
Conference on Gravitational Waves: sources and detectors, World scientific, Singapore 

Astone P., 1998, in Proceedings of the Second Edoardo Amaldi Conference, Geneva (Switerlaznd) july-1997, to appear 

Bonazzola S., Gourgoulhon E., 1996, A&A, 312, 675 

de Wild R., Nieuwkerk L.R., van Sinttruyen J.S., IEE Proceedings, 134, 659 
Frasca S., La Posta C, 1991, Nuovo Cimento, 14C, 235 

Frasca S. et al., 1992, in Proceedings of the MG6 Marcel Grossmann Meeting, World Scientific, Singapore 
Gibbons G.W., Hawking S.W., 1971, PRD, 4, 2191 

Gradshtcyn I.S., Ryzhik I.M., Tables of Integrals, Series and Products, Academic Press, San Diego 
Helstrom C.W., 1968, Statistical theory of signal detection, Pergamon Press, Oxford 



© 1998 RAS, MNRAS 000, 1-17 



All-sky search algorithms for monochromatic signals in resonant bar GW detector data 17 



Hung E.K.L., 1981, IEEE Transactions on Acoustics, Speech, and Sig Proc, ASSP-29, 163 

Hamilton W.O., 1994, in Coccia E., Pizzella G., Ronga F., cds, Proceedings of the First Edoardo Amaldi Conference, World Scientific, 
Singapore 

Kay S.M., 1990, Modern Spectral Estimation: theory and application, Prentice Hall, New Jersey 

Lobo JA., Montero M., 1997, in Miralles J. A., Morales J. A., Saez D., eds, Some Topics on General Relativity and Gravitational 

Radiation, Proceedings of the ERE-96 Conference, Editions Frontieres, Paris 
Lobo JA., Montero M., 1998, in Stela, UIB, Proceedings of the ERE-97 Conference of Palma de Mallorca, to appear 
Mauccli E., 1997, PhD thesis, Lousiana State University 

Montero M., 1997, in Buitrago J., Mediavilla E., Oscoz A., eds, General Relativity, Rclativistic Astrophysics and Cosmology, Proceedings 

of the ERE-95 Conference, World Scientific, Singapore 
Montero M., 1998, Search for monochromatic GW signals in data of the cryogenic detector Explorer, PhD Thesis, Barcelona 
Niebauer T.M. et al., 1993, PRD, 47, 3106 

Olver F.W.J., 1970, in Abramowitz M., I. A. Stegun LA., eds, Handbook of Mathematical Functions, Dover, New York 

Pallottino G.V., G. Pizzella G, 1984, Nuovo Cimcnto 7C, 155 

Papoulis A., 1984, Signal analysis, McGraw-Hill, Singapore 

Papoulis A., 1990, Probability and Statistics, Prentice Hall, New Jersey 

Roland J., Frossati G., Teyssier R., 1994, A&A, 290, 364 

Thorne K.S., 1987, in S.W. Hawking, W. Israel, eds, Three hundred years of gravitation, Cambridge University Press 
Weber J., 1969, PRL, 22, 1320 

Yip P.C.Y., 1976, IEEE Transactions on Computers, C-25, 287 

This paper has been produced using the Royal Astronomical Society/Blackwell Science LiTr^X style file. 



© 1998 RAS, MNRAS 000, 1-17 



